The present invention relates to a device for filtering an input sequence of digital signals, comprising:
(1) a first filtering stage, itself including:
(A) a first splitting circuit, provided for subdividing the input sequence into two disjoint sets of samples, comprising respectively the odd samples c0(2n+1) and the even samples c0(2n) of said sequence
(B) a first predicting circuit, provided for generating xe2x80x9cdetailxe2x80x9d coefficients d1(n)=c0(2n+1)xe2x88x92P1[ . . . , c0(2nxe2x88x922), c0(2n), c0(2n+2), . . . ], where P1 is a first linear filter and ( . . . , c0(2nxe2x88x922), c0(2n), c0(2n+2), . . . ) is the vector containing only the even samples of the input signal;
(C) a first updating circuit, provided for generating xe2x80x9capproximationxe2x80x9d coefficients c1(n)=c0(2n)+U1 [. . . , d1(nxe2x88x921), d1(n), d1(n+1), . . . ], where U1 is a second linear filter and ( . . . , d1(nxe2x88x921), d1(n), d1(n+1), . . . ) is the vector containing the xe2x80x9cdetailxe2x80x9d coefficients;
(2) at least a second filtering stage, itself including:
(D) a second splitting circuit, provided for subdividing the previously generated xe2x80x9capproximationxe2x80x9d vector into two disjoint sets of samples, similarly called c1(2n+1) and c1(2n);
(E) a second predicting circuit, provided for generating a second level of xe2x80x9cdetailxe2x80x9d coefficients d2(n)=c1(2n+1)xe2x88x92P2[ . . . , c1(2nxe2x88x922), c1(2n), c1(2n+2), . . . ], where P2 is a third linear filter and ( . . . , c1(2nxe2x88x922), c1(2n), c1(2n+2), . . . ) is the vector of xe2x80x9cdetailxe2x80x9d coefficients of even order resulting from the second splitting operation;
(F) a second updating circuit, provided for generating a second set of xe2x80x9capproximationxe2x80x9d coefficients c2(n)=c1(2n)+U2[ . . . , d2(nxe2x88x921), d2(n), d2(n+1), . . . ], where U2 is a fourth linear filter and ( . . . , d2(nxe2x88x921), d2(n), d2(n+1), . . . ) is the following vector of xe2x80x9cdetailxe2x80x9d coefficients.
This invention may be used in a wide range of signal processing applications, such as image and speech compression.
Sub-band decomposition techniques are extensively used for data coding, for instance in signal processing applications like image and speech compression. Applied for example to image compression, they allow to perform a representation of the image information as a collection of sub-images appearing at different resolutions. A popular technique for carrying out such a decomposition is the well-known wavelet transform, that allows an original input signal to be described by a set of sub-band signals. Each sub-band represents in fact the original signal at a given resolution level and in a particular frequency range.
This decomposition into uncorrelated sub-bands is implemented by means of a set of monodimensional filter banks applied first to the lines of the current image and then to the columns of the resulting filtered image. An example of such an implementation is described in xe2x80x9cDisplacements in wavelet decomposition of imagesxe2x80x9d, by S. S. Goh, Signal Processing, vol. 44, no 1, June 1995, pp.27-38. Practically two filtersxe2x80x94a low-pass one and a high-pass onexe2x80x94are used to separate low and high frequencies of the image. This operation is first carried out on the lines and followed by a sub-sampling operation, by a factor of 2. It is then carried out on the columns of the sub-sampled image, and the resulting image is also down-sampled by 2. Four images, four times smaller than the original one, are thus obtained: a low-frequency sub-image (or xe2x80x9csmoothed imagexe2x80x9d), which includes the major part of the initial content of the concerned image and therefore represents an approximation of said original image, and three high-frequency sub-images, which contain only horizontal, vertical and diagonal details of said original image. This decomposition process continues until it is clear that there is no more useful information to be derived from the last smoothed image.
In the technical report xe2x80x9cBuilding your own wavelets at homexe2x80x9d, by W. Sweldens and P. Schrxc3x6der, Industrial Mathematics Initiative, Department of Mathematics, University of Seattle, Carolina, 1995, a new way to implement a very simple example of wavelet transform is described, the so-called Haar wavelet. Consider two numbers A and B as two successive samples of a sequence (A and B have therefore some correlation), a simple linear transform allows to replace them by their average S and their difference D:
S=(A+B)/2xe2x80x83xe2x80x83(1)
D=Bxe2x88x92Axe2x80x83xe2x80x83(2)
(if A and B are highly correlated, the expected absolute value of their difference is small and can be represented with few bits). When using such a linear transform, no information is lost, since A and B can always be reconstructed:
xe2x80x83A=Sxe2x88x92D/2xe2x80x83xe2x80x83(3)
B=S+D/2xe2x80x83xe2x80x83(4)
The above-given consideration is the key for understanding the Haar wavelet transform. Consider now a signal S[n] of 2n sample values s[n, l]:
S[n]=s[n, l] 0xe2x89xa6l less than 2nxe2x88x921xe2x80x83xe2x80x83(5)
and apply the average and difference transform to each pair A=s[2l] and B=s[2l+1]. There are 2nxe2x88x921 such pairs (l=0, 1, 2, 3 . . . , 2nxe2x88x921). The results are the following:
s[nxe2x88x921, l]=(s[n, 2l]+s[n, 2l+1])/2xe2x80x83xe2x80x83(6)
d[nxe2x88x921, l]=s[n, 2l+1]xe2x88x92s[n, 2l]xe2x80x83xe2x80x83(7)
The input signal s[n], which has 2n samples, is split into two signals : s[nxe2x88x921], with 2nxe2x88x921 averages s[nxe2x88x921, l], and d[nxe2x88x921], with 2nxe2x88x921 differences d[nxe2x88x921, l], these averages and differences always allowing to recover the original signal S[n]. The averages can be seen as an approximation or a coarser representation of the original signal S[n], and the differences as the detail information needed to go from that coarser representation back to this original signal. If said original signal has some coherence (e.g. if the samples are values of a smoothly varying function), then the coarse representation closely resembles this original signal, the details being very small.
When applying iteratively the same transform to the coarser signal (by taking the average and difference), a coarser signal s[nxe2x88x922] and another difference signal d[nxe2x88x922], splitting s[nxe2x88x921], are obtained. When repeating this operation n times, a Haar transform is implemented (that can be thought of as applying a Nxc3x97N matrix (N=2n) to the signal S[n]).
A new way of looking at the Haar transform, also called the xe2x80x9cladderxe2x80x9d of lifting scheme, can then be presented, the novelty lying in the way the difference and the average of two numbers A and B can be computed. If one wants to compute the whole transform in-placexe2x80x94i.e. without using auxiliary memory locations, by overwriting the locations that hold A and B with the values of respectively S and D-, it cannot be done with the previous formulas (1) and (2): when using them and storing S and D in the same location as A and B respectively, it would lead to a wrong result.
Another implementation in two steps is preferred:
first the difference D=Bxe2x88x92A is computed and stored in the location for B (which is therefore lost, thus preventing from computing directly S=A+B);
B being lost, the average is then found by using A and the newly computed difference D, according to the formula S=A+D/2 which gives the same result, and this result is stored in the location for A.
The advantage of such a splitting into two steps is that it is possible to overwrite B with D and A with S without requiring any auxiliary storage. This particular scheme can be described in more detail. Consider a signal S[j ] with 2j samples, which has to be transformed into a coarser (or xe2x80x9capproximationxe2x80x9d) signal S[jxe2x88x921] and a xe2x80x9cdetailxe2x80x9d signal D[jxe2x88x921]. A typical case of a wavelet transform built through the above scheme (see the technical report previously cited) consists of three steps, each of which is described in more detail, in relation with FIG. 1:
a splitting step SPL: the signal S[j] is divided, in a splitting stage 11, into two disjoint subsets of samples, one group consisting of the even indexed samples s[2l] and the other one consisting of the odd indexed samples s[2l+1] (this splitting into even and odd samples is called the Lazy wavelet transform). The following operator can therefore be built:
(even [jxe2x88x921], odd [jxe2x88x921])=split[s[j]]xe2x80x83xe2x80x83(8)
a predicting step PRED: if the signal S[j ] has a local correlation structure, these even and odd subsets are highly correlated, which allows to predict, in a predicting stage 12, one of the two sets with a reasonable accuracy if the other one is given (in the present case, the even set is always used to predict the odd one: an odd sample s[j, 2l+1] uses its left neighbouring even sample s[j, 2l] as its predictor, for determining the difference d[jxe2x88x921, l]=s[j, 2l+1]xe2x88x92s[j, 2l] between the odd sample and its prediction). More gene other subband transforms can be implemented using a prediction of the odd samples based only on the previous even ones. Denoting by P the prediction operator, the detail signal is then obtained as the difference between the actual odd sample, at the output of the splitting stage 11, and its prediction, at the output of the predicting stage 12, and can be written:
d[jxe2x88x921]=odd[jxe2x88x921]xe2x88x92P[even[jxe2x88x921]]xe2x80x83xe2x80x83(9)
where xe2x80x9cevenxe2x80x9d represents all even samples previous to jxe2x88x921.
an updating step UPD: one of the key properties of the coarser signals is that they have the same average value as the original signal, i.e. the quantity                     S        =                              2                          -              j                                ⁢                                    ∑                              l                =                0                                            2                                  j                  -                  1                                                      ⁢                          xe2x80x83                        ⁢                          s              ⁡                              [                                  j                  ,                  l                                ]                                                                        (        10        )            
is independent of j, which has for result that the last coefficient s[0,0] is the DC component or overall average of the signal. The updating step, carried out in an updating stage 13, ensures this by letting:
xe2x80x83s[jxe2x88x921, l]=s[j, 2l]+d[jxe2x88x921, l]/2xe2x80x83xe2x80x83(11)
Substituting this definition, it can easily be verified that:                                           ∑                          l              =              0                                      2                              j                -                1                                              ⁢                      xe2x80x83                    ⁢                      s            ⁡                          [                                                j                  -                  1                                ,                l                            ]                                      =                                            ∑                              l                =                0                                            2                                  j                  -                  1                                                      ⁢                          xe2x80x83                        ⁢                          (                                                s                  ⁡                                      [                                          j                      ,                                              2                        ⁢                        l                                                              ]                                                  +                                                      d                    ⁡                                          [                                                                        j                          -                          1                                                ,                        l                                            ]                                                        ⁢                                      /                                    ⁢                  2                                            )                                =                                    1              ⁢                              /                            ⁢              2              ⁢                                                ∑                                      l                    =                    0                                                        2                                          j                      -                      1                                                                      ⁢                                  xe2x80x83                                ⁢                                  (                                                            s                      ⁡                                              [                                                  j                          ,                                                      2                            ⁢                            l                                                                          ]                                                              +                                          s                      ⁡                                              [                                                  j                          ,                                                                                    2                              ⁢                              l                                                        +                            1                                                                          ]                                                                              )                                                      =                          1              ⁢                              /                            ⁢              2              ⁢                                                ∑                                      l                    =                    0                                                        2                    j                                                  ⁢                                  xe2x80x83                                ⁢                                  s                  ⁡                                      [                                          j                      ,                      l                                        ]                                                                                                          (        12        )            
which allows to define an operator U of the form:
s[jxe2x88x921]=even[jxe2x88x921]+U[d[jxe2x88x921]]xe2x80x83xe2x80x83(13)
In the case of the Haar transform, U[d[jxe2x88x921]] is simply (d[jxe2x88x921])/2. More complicated updating operators exist for other subband transforms.
The subtractor and the adder are designated by the references 14 and 15. All these computations are made in-place: the even locations can be overwritten with the averages and the odd ones with the details. A C-like implementation of these computations may be given by:
(odd[jxe2x88x921], even[jxe2x88x921]):=Split (s[j])
odd[jxe2x88x921]xe2x88x92=P(even[jxe2x88x921])
even[jxe2x88x921]+=U(odd[jxe2x88x921]
It has then been shown that any linear wavelet transform can be implemented with such kind of scheme, where the scale transition is achieved by applying a second predict/update stage to s and d signals, thus constituting a filter bank. An important characteristic of the filter bank which is employed should be the maximization of the energy compaction of the input, which drastically influences the final coding gain. In this context, efforts have been made on the design of filter banks that are optimally adapted to the input signal statistics.
These filter banks are generally non-varying decompositions, which act in the same way whatever the input image. It then appeared as possible to design adaptive filters, changing with the varying characteristics of the input, but, till now, only the first step of these schemes was adapted to changing input, such as described for instance in xe2x80x9cLinear/nonlinear adaptive polyphase subband decomposition structures for image compressionxe2x80x9d, by xc3x6.N. Gerek and A. E. Cetin, Proceedings of the 1998 IEEEE International Conference on Acoustics, Speech and Signal Processing, vol.III, May 12-15, 1998, Seattle, Wash. (USA), pp. 1345-1348.
It is therefore an object of the invention to propose, with respect to such an implementation, an optimized adaptive filtering device, allowing a complete adaptation to varying characteristics of the input.
To this end, the invention relates to a device such as defined in the introductive paragraph of the description and which is moreover characterized in that the first updating circuit U1 and the second predicting circuit P2 are associated within an iterative cross-optimization operation carried out by using the second filter U1 for optimizing the third filter P2 and then said third filter P2 for optimizing said second filter U1.
The advantage of this structure, adapted to the non-stationarities in the input signal, is that it preserves the perfect reconstruction property, without respect of the filters employed in each concerned pair of filtering stages. This is a major advantage in image coding schemes, because lossless coding is achieved.
Preferably, the iterative optimization operation aims at minimizing the variance of the detail coefficients at the output of each decomposition level. The efficiency of the encoding operation of each detail sub-band is indeed increased if its variance has been minimized. At the last level, the optimization aims at minimizing the variance of the approximation coefficients.